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We investigate the influence of tfie electron-pfionon coupling in the one-dimensional spinless Hol- 
stein model at half-filling using both a recently developed projector-based renormalization method 
(PRM) and an refined exact diagonalization technique in combination with the kernel polynomial 
method. At finite phonon frequencies the system shows a metal-insulator transition accompanied 
by the appearance of a Peierls distorted state at a finite critical electron-phonon coupling. We 
analyze the opening of a gap in terms of the (inverse) photoemission spectral functions which are 
evaluated in both approaches. Moreover, the PRM approach reveals the softening of a phonon at 
the Brillouin-zonc boundary which can be understood as precursor effect of the gap formation. 

PACS numbers: 71.10.Fd, 71.30.-fh 

I. INTRODUCTION 

Despite the many years of study of the electron-phonon interaction in metallic systems, there remain fundamental 
problems that yet have to be resolved. Especially systems which suffer strong electron-phonon coupling in conjunction 
with strong electron-electron interaction are in the center of present interest. Examples are cuprate high-temperature 
superconductorsii^, colossal magnetoresistive manganites^, or metallic alkaline-doped Ceo-based compounds. Further- 
more, in a wide range of quasi-one-dimensional materials, such as MX chains, conjugated polymers or organic charge 
transfer complexes*, the itinerancy of the electrons strongly competes with the electron-phonon coupling which tends 
to establish e.g. charge-density-wave structures. Then, in particular at half-filling, Peierls insulating phases may be 
energetically favored over the metallic state. Many interesting questions arise not only with a view to the associated 
metal to insulator transition but also concerning the form of the single-particle excitation spectra well below and 
above the transition. At present, there is a clear need of reliable theoretical methods to tackle these problems in 
terms of minimal microscopic models. 

In this paper we study as the perhaps simplest realization of a strongly coupled electron-phonon system, the so- 
called spinless Holstein model, describing the local interaction between dispersionless longitudinal optical phonons 
and the density of the electrons at a given lattice site = cjci: 

n = -i^(4c, +h.c.)+^o5Z blh+gJ2ibl+h)n,. (1) 

Here, {bl) denote fermionic (bosonic) creation operators of electrons (phonons), i is the Wannier site index. The 
electron-phonon coupling constant and frequency of the Einstein mode is given by g and loq, respectively. Of physical 
concern for most applications are relatively small values of the adiabaticity ratio a = tuo/t <C 1, even though the 
anti-adiabatic limit a ^ 1 is also a useful point of reference for an overall understanding of the physics of the Holstein 
model. 

It is well-known that the one-dimensional Holstein model of spinless fermions at half-filling has a quantum phase 
transition from a Luttinger liquid (metallic phase) to an insulating phase with charge-density wave long-range 
orde r^i^i^ . In the past, a large number of different analytical and numerical methods have been applied to the Holstein 
model, in particular to determine the phase boundary between metallic and insulating behavior for the half-filled band 
case. Much of the work is restricted to the one-dimensional case. Mainly ground state properties were investigated 
by means of strong-coupling expansions^, variational^ and renormalization group&S approaches, as well as world-line 
quantum Monte Carlc^ and Green's-function Monte Carlc^" simulations. More recently, exact diagonalization^^ (ED) 
and density matrix renormalization group2ii2ii^ii4 techniques were applied. The metal-insulator transition is accom- 
panied by the appearance of a gap in the electronic spectrum which, however, can best be observed in the k-dependent 
one-particle spectral functions. In a recent dynamical mean field treatmen^i^ in conjunction with a numerical renor- 
malization group approach it was suggested that the opening of the electronic gap is accompanied by the appearance 
of a low-energy phonon peak in the total phononic spectral function. In fact, we shall show in this paper that the 
phonon modes at the Brillouin-zone boundary become soft. This can be understood as a precursor effect of a lattice 
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instability leading to a Peierls state for electron-phonon coupling strength g larger than some finite critical value gc- 
Moreover, we shall evaluate electronic one-particle spectral functions which should show the opening of the gap at 
the quantum phase transition. The evaluation of Luttinger parameters is beyond the scope of interest of this paper. 

One of the main aims of the paper is to show that a newly developed projector-based renormalization method 
(PRM) for many-particle Hamiltonians^^ can be applied to the spinless Holstein model, though an extension to the 
case with spin would also be possible. In principle this method is applicable to infinitely large systems, to the whole 
parameter regime of the electron-phonon coupling g/too, to any finite filling, and to any finite spatial dimension. 
In the present paper we restrict ourselves to the case of one dimension and to half-filling in order to study the 
transition from the metallic to the insulating phase. It will turn out the Peierls instability associated with the metal- 
insulator transition can very well be described in the present formalism. Since the application of the PRM, however, 
is accompanied by some approximations, the reliability of the technique will be tested by comparing its results with 
unbiased data from exact diagonalization (ED). In particular, the k-dependent one-particle spectral functions from 
(inverse) photoemission will be analyzed within both approaches. As is seen, the formation of a charge gap at some 
critical value gc of the electron-phonon coupling is related to the softening of phonon modes at the Brillouin-zone 
boundary, at least in the adiabatic regime. From the ED data we get additional valuable insights into the behavior 
of the wave-vector resolved spectral function. Moreover, a detailed characterization of the ground state, e.g. by the 
phonon distribution is possible. We note that the present approach is not restricted to large dimensions as is the case 
of a recent DMFT approacbi^. 

The paper is organized as follows. In Sect. II we outline the PRM approach^^. The renormalization equations for 
the model parameters are derived non-perturbatively and general expressions for the one-particle spectral functions 
are discussed. In Sect. Ill we discuss our exact diagonalization technique for the calculation of single-particle spectral 
functions of coupled electron-phonon systems. Sect. IV presents the findings obtained within both approaches bf for 
the (I)PE spectra and phonon renormalization, in particular with respect to the metal to Peierls insulator transition 
in the ID Holstein model. The main results will be summarized in Sect. V. 



II. RENORMALIZATION OF THE HOLSTEIN MODEL 



A. Projector-based renormalization method 



The PRM^^ starts from a decomposition of a given many-particle Hamiltonian Ti, into an unperturbed part Ti.o and 
a perturbation Tii 



where we assume that the eigenvalue problem of Ho is solved 



(2) 
(3) 



The decomposition of Ti. into Tio and Hi should be done in such a way that Hi contains no part which commutes 
with Ho- Thus, Hi gives rise to transitions between eigenstates of Ho with different eigenenergies. The presence of 
Hi usually prevents the exact solution of the eigenvalue problem of the full Hamiltonian. Let us define a projection 
operator Pa by 



El 



A 



m(°)\e(A-|i?(°)-£;(°)|). 



(4) 



Pa is a super-operator which acts on ordinary operators A of the unitary space. It projects on those parts of A 



p(o) 



(0) 



which are formed by transition operators (to'"'! with energy differences 

A, where A is smaller than the cutoff A of the original model. Note that in Eq. Q neither 
be low-energy eigenstates of Ho- However, their energy difference has to be restricted to values < A. Furthermore we 
define the projector 



less than a given cutoff 
''■°'') nor Ito^"-*) have to 



Q 



A — 



(5) 



on the high-energy transitions larger than the cutoff A, (QaPa = 0). 

Now we want to transform the initial Hamiltonian H of Eq. into an effective Hamiltonian H\ which has no 
matrix elements belonging to transitions larger than A. This will be achieved by an unitary transformation 



Hx 



He 



(6) 
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where the generator X\ of the transformation has to be anti-Hermitian, i.e. X\ = —X\. Note that the effective 
Hamiltonian T-Lx has the same eigenspectrum as the original Hamiltonian 7i. The generator Xx has to be chosen such 
that Tix has no matrix elements between states belonging to transitions larger than A. Hence the condition 

QaWa - (7) 

has to be fulfilled. Eq. Q will be used below to specify X\. 

Instead of eliminating all high-energy excitations in one step a sequence of stepwise transformations will be used 
in the following. Thus, in an infinitesimal formulation, the projector-based rcnormalization approach yields renor- 
malization equations for the parameters of the Hamiltonian as function of the cutoff A. In that aspect the present 
approach resembles Wegner's flow equation method^'' and the similarity transformation introduced by Glatzek and 
WilsoniS.. To find the rcnormalization equations we proceed as follows. We start from the renormalized Hamiltonian 

'Hx = Ho, A + Hi, A , (8) 

where in Tlx all excitations with energy differences larger than A have been eliminated. Next we further integrate out 
all excitations inside an energy shell between A and a smaller cutoff (A — AA) where AA > 0. The new Hamiltonian 
'H(A-AA) is given by 

H(A-AA) = e^^'^^ Hx e-^^'^^ , (9) 
where Xx,a\ is the generator for the transformation from A to (A — AA). Similar to Q it has to fulfill the condition 

Q(A-AA) H(A-AA) = 0. (10) 

Note that there are two strategies to exploit Eq. 1)10(1 in order to determine the generator Xx,a\ of the unitary 
transformation (EJ. The most straightforward route is to analyze Eqs.Q and in perturbation theory as it was 
done in Rcfs. 16 and 19. Here, we want to perform the rcnormalization step from A to (A — AA) in a non-perturbative 
way which recently has been applied to the periodic Anderson model by two of the authorsSS. 

Eqs. (PJ and ((10(1 can be used to derive difference equations for the A dependence of the parameters of the Hamilto- 
nian. They will be called rcnormalization equations. The solution depends on the initial values of the parameters of 
the Hamiltonian and fixes the final Hamiltonian Tlx in the limit A — > 0. Note that the final Hamiltonian only consists 
of an renormalized unperturbed part ?io,(A^o)- The interaction 7ii.(A^o) completely vanishes since it was used up in 
the rcnormalization procedure. 



B. Application to the Holstein model 

Let us start by formally writing down the effective Hamiltonian Tix = e^^Tie^-^^ for the spinless Holstein model 
after all excitations with energy differences larger than A have been eliminated, 

Ha — Ho, A + Hi, A 

with 

Ho, A = £k,A cj^Ck + CJq,A b\h^ + Ex , (11) 

k q 

Hi,A ^ 5: ^ Pa [h\cic^.+^ + 6,c^,Ck) , (12) 



where Fourier transformed operators 



k,q 



1 3 



were introduced. Due to the rcnormalization processes the one-particle energies ek,A and the phonon frequencies o^q.A 
in l(ll() now depend on the cutoff A. Moreover, the phonon energies acquire a dispersion due to an effective interaction 
of lattice vibrations at different sites via the coupling to electronic degrees of freedom. Also the electron-phonon 
coupling constant 5k. q. A now depends on wave vectors k, q and the cutoff A. Ex is an additional energy shift. Finally, 
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the projector in Eq. (|12|l guarantees that only those excitations survive in Hi^x which have energies (with respect 
to Ti-Q x) smaher than A: 

Pa (&^4ck+q + Vl+qCk) = (13) 

= 9 (A - |a;q,A + £k,A - ek+q,A|) (^qCk^k+q + ^qCfc+qCk) ■ 

The initial values of the original model (with cutoff A = A) are 

£k,(A=A) = ^k, 'i^q,(A=A) = ^0 , 5k,q,(A=A) = 5 J i?(A=A) = ■ (14) 

In the next step we determine the renormalized Hamiltonian 7i(A-AA) — e^^'^'''H\e~'^^ '^^ by eliminating all 
excitations within an additional small energy shell between (A — AA) and A. For the explicit form of the generator 
^A,AA of the unitary transformation we make the following ansatz 

^A,AA - ^^Bk.q,Aek,q(A,AA) (fo^4ck+q - Vl+qCk) , (15) 



k,q 

where 0k,q(A, AA) is the product of two 0-functions 

ek,q(A,AA) = e(A- K,A+ek.A-£k+q,A|) (16) 

X0 [kq,(A-AA) + £k,(A-AA) - £k+q,(A-AA)l - - ^^)] ■ 

The operator form of Xa,aa is suggested by its first order expression which is easily obtained by expanding (jSJ in 
powers of Tii and using ((117)) . (cf. Ref. H^ . The yet unknown prefactors -Bk,q,A will be specified later and depend 
on A. It will turn out that i?k.q.A contains contributions in all powers of the electron-phonon coupling g. The two 
0-functions in Hlt)|) confine the allowed excitations to the energy shell AA. 

The coefficients i?k.q.A will be fixed by the condition H10() . First, we have to carry out the unitary transformation 
© explicitly: 

W(A-AA) = ^k.A e^^.^vj^Cke-^^-^^ + ^ ujq,x e^^.^^6]jfeqe-^^.^^ + Ex + 

k q 
k,q ^ 

The transformations for the various operators in (|17|l has to be done separately. For instance, the transformation for 
cj^Ck reads 

= E I f ^'""l"^^'^^^ (cOs(Sk-q,q,Ay2^) - l) {(1 - "£-q,A + <.a) ^Ck 
q ^ L ^ 

-("k.A + '^q.A) Ck-qCk-q} 

ek-q,q(A,AA)(?l^ ;, J J 

H Z ■ { C0S(Bk-q,q.AV2"-q,A) ^ 1 > o' 6q 

n-q,A ■ J 

+ ^^^^^#^Bin(Bk-q,q,Av/2^) fecLqCk + h., 

where we have defined 



(4ck)o , = -T ' <,A = (&q^>0,A = _ 1 ' (19) 



'k.A \-k-K/n.A g;3£k,A + 1 ' '"q.^ \"q"q/o,A g/3c^q,A _ 1 

and 



nq,A = l + 2<,;, (20) 
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(for details see Appendix^ . In deriving (|18|l an additional factorization approximation has been used in order to keep 
only operators of the structure of those of Eqs. Ullfl and H12I) . In principle, the expectation values are best defined 
with the equilibrium distribution of Tlx since the renormalization step was done from T^a to 7i(A-AA)- However, for 
simplicity we shall evaluate these quantities with the unperturbed Hamiltonian 7io,A- Similar expressions to (|18|) are 
also found for the transformations of the remaining operators of (|17|l . 

One should note that influence of the above factorization approximation together with the choice of the equilibrium 
distribution is not well controlled as long as the influence of the additional fluctuation terms are not considered. 
However, we believe that the renormalization of Tlx as it is used here leads to a proper description of the influence of the 
electron-phonon interaction in the Holstein model. The fluctuation terms which appear as additional renormalization 
contributions in the approach could give rise to additional terms both in 7io,A and Ti-i^x- The renormalization of the 
new coupling parameters would have to be investigated in the present renormalization scheme as well. 

In the next step we determine the parameters -Bk,q,A- For that purpose we use the condition (|10|) . By inserting 
^(A-AA) from H17|l into H10() and by use of IjlSI) and the remaining transformations, we find 



ek,q(A,AA) Bk.q.A = ek,q(A,AA) 



: arctan 



I2h, 



q,A 



5k,q,A 



N 



(21) 



Let us point out that -Bk,q,A is determined by (|21|) only for the case that the excitation energies LOq^x + £k,A ~ ^k+q.A 
fit into the energy shell given by 6k,q(A, AA). For all other excitations -Bk,q,A can be set equal to zero. Therefore we 
shall use the following expression for -Bk.q.A^ 



1 



k,q,A 



: arctan 




for ek,q(A, AA) = 1 
for ek,q(A, AA) = . 



C. Renormalization equations for the Holstein model 



Next, we derive the renormalization equations for the parameters of the Hamiltonian. For that purpose we compare 
the renormalization ansatz for Ti^x-AX) i 

W(A-AA) = X^^'kXA-AA) 4ck + X!^q.(A-AA) feq&q + £^(A-AA) (22) 
k, q 

+ E P(^-AA) (fel4ck+q + 6qcUqCk) 

k,q ^ 

[see Eqs. and H12|) ]. with the expression that is obtained from (|17|l after (|18l) and the corresponding transformations 
have been inserted. Comparing in both equations the coefhcients of the operators c|j.Ck, &q&q, and (foj^cj^Ck+q + h.c), 
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we find the following relations between the parameters at cutoff A and (A — AA): 

£k,(A-AA) - £k,A = (23) 
^ X! ^^"^'^ ^ '^^^ |cOs(Bk,q,AV2"-q,A) - 1} ('^q,A + £k,A - £k+q,A) ©k,q(A, AA) 



q "1^^ 



\l ("k+q.A + "q,A) sin(Bk,q,A\/2nq,A)5k,q,Aek,q(A, AA) 



q 



!i_q^ ELA I COs(-Bk-q,q, A a/ 2nq, a) - l| (Wq,A + Ek-q.A " ek,A) 6k-q,q(A, AA) 



^ \l Nn ^ "k-q,A + "q,A) Sin(Bk-q,q,AA/2nq,A)ffk-q,q,A6k-q,q(A, AA) , 



q V '^'^ 



t^q.CA-AA) - Wq,A = (24) 
"q,A 



^ k+q^A k^ |cOs(Sk,q,A V2nq,A) - l| (t^q.A + ^k.A - ek+q,A) 6k,q(A, AA) 



^ \l Nn ("''+1^^ ~ "^,a) sin(Bk,q,A\/2?Jq,A) 5k,q,Aek,q(A, AA) , 



k 



q,A 



.9k,q,(A-AA) ~ 5k,q,A — (25) 



N 



2h. 



q,A 



(tjq,A + ek,A - £k+q,A) sin(i?k,q,A V^^q^) 0k,q(A, AA) 



-5k,q,A |cOs(Bk,q,Av/2?^q,A) ^ l| 0k,q(A, AA) . 



A corresponding expression can also be found for the renormalization of the energy shift £'(a-aa)- The above Eqs. 123(1 
to H25|l describe the renormalization of the parameters if the cutoff is reduced from A to (A — AA). 

Note that the expression ((21|l for i?k,q,A and thus the renormalization equations are non-pcrturbative in (7k,q,A and 
are not restricted to some low order. However, a prefactor 1/VN enters the expression (|21ll for i?k.q,A due to the 
factorization approximation discussed in Appendix A. Thus, in the thermodynamic limit N 00, i?k,q,A becomes 
linear in g and the renormalization relations (|23|l to (|25|l become quadratic in the electron-phonon coupling. In the 
numerical evaluation, however, when a fixed number N is taken, the excitation energy (wq,A + £k,A — £^k+q,A) will 
become very small for reducing the cutoff A — > and thus the expansion of -Bk,q,A to linear order in g may break 
down. In that case the full equations pifl and l(23|l to l(25|) have to be taken. 

The overall renormalization starts from the cutoff A = A of the original model and proceeds down to zero cutoff 
A = 0. At A = the completely renormalized Hamiltonian Ti. := 7i(A^o) becomes an effectively free model and reads 

^ = ^ek4ck + ^tiq6lj6q + ^, (26) 

k, q 

where we have defined Ck — ek.(A^o)j<^q = ^q.(A^o)j f-i^d E = _E(a^o)- Foi' A — > the electron-phonon coupling Ti.i^\ 
has vanished due to the 0-function in p3|l . All excitations of the electron-phonon interaction were used to renormalize 
the parameters of the Hamiltonian. 

The results from the numerical evaluation of the renormalization equations will be given in Sect. IV. In particular 
we then compare the spectral functions calculated within the PRM approach and by exact diagonalization. The 
electronic one-particle spectral functions A^{uj) and A^{uj) are defined by 

AU^) - 2^ I (ck(t)4)e-*dt, A^{uj)^—I {ci{t)c^)e^-'dt. (27) 



— 00 -J —00 
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The function A^{lu) describes the creation of an electron k at time zero and the annihilation at time t whereas in 
A^{u!) first an electron is annihilated and at time t the electron is created. As is well-known these quantities can 
be measured by inverse photoemission (IPE) and by photoemission (PE). To evaluate A^{uj) and A^{uj) within the 
PRM approach it is necessary to unitary transform not only the Hamiltonian but also the operators Ck and cj^. This 
follows from the fact that the trace of any operator quantity is invariant under a unitary transformation. Thus we 
have to evaluate 

Ki-^) = ^ / ((ck)AW(4)A),e-*dt, (28) 

where the expectation value and the time dependence are defined with respect to the A-dependent Hamiltonian Hx- 
Moreover, {cl)\ = e-^^cl^e"-^^ . A similar expression to Eq. (|28|l is also valid for A^{uj). 

In analogy to the evaluation procedure for the renormalization equations of Tix we make the following ansatz for 
the A-dependence of the operators (c|^)a and (ck)A, 

(4)a - ak,A Cl+J2 (/5k.q.A 4+q6q + 7k,q,A 4.^6^) , (29) 

1 

(ck)A = [(ck)!]^; with A-dependent parameters akA: /3k,q,Aj and 7k, q, a- The parameter values for the original model 
are 

ak,(A=A) = 1 , /3k,q,(A=A) = , 7k,q,(A=A) = . (30) 

Due to H29|l the following sum rule must hold 

1 = |ak,A|' + ^{|/?k,q,A|'K+q,A+<A) + 

q 

+ |7k,q,Ap(l+<A-<-q.A)} , (31) 

which follows from the commutator relations. Note that a factorization approximation was used on the right hand side 
of H31|) . In analogy to the former approach for £k.A, Wq.A, and gk.q.A, we find the following renormalization equations 
for the parameters akA, /^k,q,A and 7k, q, a: 

ak,A--AA - ak,A = (32) 

^ E {^°^ ( V"k+q,A + <,A ^k,q,A) - 1} «k,Aek,q(A,AA) 

q 

+ E (V^l - "k-q,A +<A 5k-q,q,A) - l} «k, A Ok-q.q ( A, A A) 

q 

+ E \/"k+q,A+<A Sin (y^n^+q,A+<A 5k,q,A) /3k,q, A ek,q ( A, AA) 

q 

- E \/l-"k-q,A+<A sin (\/l-"£-q,A + <A ^k-q,q,A) 7k,q, Aek-q,q( A, AA) , 

q 



and 



/3k,q,A-AA — /3k, q, A = 
1 



(33) 



"k+q,A + '^q.A 



sin ( A/?lk+q,A + "q,A ^k,q,A j ak,Aek,q(A, AA) 



+ 



I cos 



Sk,q,A) /3k,q,Aek,q(A, AA) , 



7k,q,A-AA — 7k, q, A — 
1 



(34) 



1 - nl 



: sin ( A /I - ?^k_q,A + "q,A -Sk-q,q,A j ak,A6k-q,q(A, AA) 



k-q,A I "'q,A 



I cos ^ 



1 - "k-q,A + 



Bk-d.d.x) - 1} 7k,q,Aek-q,q(A,AA), 
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The renormalization equations (|31|l and H33|l can again be integrated numerically by reducing the cutoff A stepwise 
down to A ^ 0. Together with the Eqs. (|23|l to H25I) the spectral functions (|27|) can be determined. 

q 

+ 7k,q - '^q - Ek-q) (1 - + U^) , (35) 

q 

^ki^) = "k '5(w-£k)r~lk+X!^k,q^('^q-£k+q+'^) r'lk+qll+^Jq) 

q 

+ 7k,q '5(^ - tl'q - £k-q) K-^K ■ (36) 

q 

As before, the quantities with tilde denote the parameters at cutoff A 

&k = ak,(A^O) /3k,q = /3k,q,(A^0) 7k,q = 7k,q,(A^0) • (37) 

In deriving (|35() and H36() it was exploited that the completely renormalizcd Hamiltonian H is diagonal [cf. Eq. (|26l) ]. 
Note that the first parts in A^{llj) and (lu) describe the coherent one-electron excitations which correspond to those 
of a free electron gas. The two remaining contributions are incoherent excitations due to the coupling to phonons 
and are given by a q-suni over excitations £k±q T <^q- Finally we give a sum rule which follows from the frequency 
integral over the sum of the two spectral functions 

1 = r du{A+{u)+A^{u:)) = 

J — OO 

= «k + E /3k,q("k+q + O + E 7k,q(l + ^ " "k-q) ■ (38) 

q q 

Note that Eq. H38|l is equivalent to the former relation H31() . The explicit evaluation of the spectral functions will be 
performed in Sect. IV below. 



III. EXACT DIAGONALIZATION OF ELECTRON-PHONON MODELS 



In principle exact diagonalization techniques allow the analysis of ground-state and spectral properties of microscopic 
models free of any approximations. However, the vast Hilbert space dimensions restrict ED studies to rather small 
system sizes, even if state of the art supercomputers are used. Usually the dimension of the matrix involved in 
the diagonalization is reduced by exploiting lattice or spin symmetries'^. Unfortunately, for the Holstcin model the 
Hilbert space associated to the phonons is infinite even for finite systems. Thus a well controlled truncation procedure 
has to be developed^LSS. In our approach the maximum number of phonons per state (not per lattice site!) is fixed 
(M). Then the dimension of the phononic subspace is Dph — ^^^^ffr'- 

For the Holstein model computational requirements can be further reduced. It is possible to separate the symmetric 
phonon mode, Bq = bi, and to calculate its contribution to H analytically. For the sake of simplicity, we restrict 

ourselves to the ID case in what follows. Using the momentum space representation of the phonon operators and 
introducing the polaron binding energy Sp — g'^/ujo the original Holstein Hamiltonian reads 

n = -tJ2icUj + h.c.) - E(^-Q. + Bq>Q. + E K.^Q. (^9) 

with 

= Uj,i&I, SQ^=U],A=U_j,A, (40) 

i 

where Uj i = exp {iQjRi}. Qj (i?^) denote the allowed momentum (translation) vectors of the lattice. The 
phononic Q = mode couples to no = ^= which is a constant if working in a subspace with fixed number of 
electrons. Thus the Hamiltonian decomposes into H — Ti. + 'Hq=q, with 

nQ=o^~^^iBl+B„)no+LOoBlB„. (42) 
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Since [H , 'Hq=o] = holds the eigenspectrum of Ti, can be built up by the analytic solution for 'Hq=o and the 
numerical results for Ti.. Using the unitary transformation 

S(iVeO = cxp \-^MiBl -bA, (43) 

and introducing a shift of the phonon operators {Bq Bq + ^^y^^), we easily find the diagonal form of T-Iq=o 

nQ=o - ^^bIb, ~ . (44) 

It represents a harmonic oscillator with eigenvalues and eigenvectors: 

Ej = u^l-ej,^, (45) 
10 = ^(4)'|0). (46) 

The corresponding eigenenergies and eigenvectors of 'Hq=o are Ei = Ej and 

\l,N,i)^S\N,i)\l), (47) 

respectively. That is, in the eigenstates of the Holstein model a homogeneous lattice distortion occurs. Note that the 
homogeneous lattice distortions are different in subspaces with different electron number. Thus excitations due to 
lattice relaxation processes show up in the one-particle spectral function. Finally, eigenvectors and eigenenergies of 
7i can be constructed by combining the above analytical result with the numerical determined eigensystem (En^"'^', 
|^r='))) ofTY: 

EiN.O = Eir^'^+.ol-e,^, (48) 

iC"^) = I'^r='^>®l^,^e0- (49) 

The direct product notation for the eigenvectors of Ti implies that also the one-particle spectral functions, e.g. the 
single-particle spectral function 

71,1 

can also be decomposed in a numerical and an analytical contribution: 

A+{uj) = J2 + 1, Nei) A+{cj - [iuol - ep{2Nei + l)/iV]) . (51) 

I 

Using ED in combination with kernel polynomial moment expansion and maximum entropy methods22i22iS4, we first 
compute numerically the one-particle spectral function for 7Y, 

i+ (ci) = E K'^^^^'^'^l4l'^^'^>P - iM!"^'^'^ - ^r^'^]) , (52) 

n 

which does not include the effects of the Q — phonon mode. In a second step, the final spectrum is constructed 
from Eq. I|51|l by shifting A'^{Cj) by multiples of the bare phonon frequency and using the weight factort 
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p{l',Ki,N,i) = |(7V^,,r|0,iVeOP = ^ e-^ (53) 

where x = — N^i)'^. The p{l' , N^i, Nei) quantify the relaxation process of the lattice distortion discussed 

above. Thus, each excitation in A'^{Cj) splits up into a band of peaks separated by ujq which contains the total weight 
of the original excitation [J2i p{l,^ei + = !]■ 
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Fig-fflshows the contribution of the Q = phonon mode to ^^(w). The Q = excitations are separated by cjq and 
interfere at higher energies with excitations forming broad bands. In Fig^ e.g., five repHcations of the lowest peak 
are visible (with decreasing height and weight). The integrated spectral density reflects the redistribution of spectral 
weight to higher energies but conserves the spectral weight when integrating the Q ^ peaks together with the 
corresponding Q — side bands. Note that this behavior is observed in a finite system only, i.e., p{l' > 1, N^i ± 1, N^i) 
vanishes in the thermodynamic limit N oo. 

In our finite cluster diagonalization, however, the separation of the phononic Q = mode significantly reduces 
the computational requirements. First, the dimension of the matrix to be diagonalized is reduced by a factor of 
1 + ^ because only — 1 instead of N independent modes have to be considered. Second, the Q mode 

takes into account at least (BqBq) — ^-j^ phonons in the ground state already. Thus, the maximum number 

of phonons allowed per state, M, can be chosen much smaller for H than for the full Hamiltonian H in order 
to achieve the same level of convergence. The relevance of the separation of the phononic Q = mode for our 
numerical work becomes evident if one compares the dimension of the phononic Hilbert space used for the calculation 
in Fig. U] {Nei = 4, iV - 1 = 7, Af = 28, Ep/wo = 6) Dph = 6.7 x 10^ with that required for the full Hamiltonian 
{N — 8, M — 28 + 12) Dpii — 3.77 x 10^, which is about two orders of magnitude larger. 

IV. PEIERLS TRANSITION IN THE SPINLESS FERMION HOLSTEIN MODEL 

In this section we apply the PRM and ED techniques outlined in Sees. II and III to the investigation of the metal 
insulator transition in the Holstein model. While both approaches are valid for any dimension in principle, we restrict 
ourselves to the one-dimensional case and spinless fermions, mainly because of the serious memory restrictions within 
the ED calculations. On the other hand, focusing on the half-filled band case N^i = N/2, we know that the Peierls 
instability, we are interested in, is most pronounced in low-dimensional systems. Throughout the numerical work 
we use the bare tight-binding electron dispersion Sk — ~2t cos ka (lattice constant a ~ I), and an Einstein phonon 
frequency ujo/t — 0.1. The Fermi level is defined to be the energy zero level and the temperature is T = 0. In the 
following we shall vary the dimensionless electron-coupling constant g/t or, equivalently, Sp/t [where due to g — ^£pUjQ 
we have Sp/t = [g /t)'^{t/LOQ)]. 

A. ED results 

In order to investigate lattice dynamical effects on the Peierls transition we first analyze the spectral density of 
single-particle excitations associated with the injection of an electron with wave number K, A^{u;) (IPE), given by 
Eq. l(5n|) . and the corresponding quantity for the emission of an electron, A'^{lo) (PE). In A~^{uj) the destruction 
operator ck connects the ground state of the Holstein model with N^i electrons (I^Aq^"''')) to all excited states of the 

system with N^i — 1 electrons dV'!^"' ^^))- Of course, the exact determination of A^{ijj) for coupled electron- phonon 
systems is an tremendous numerical task, which requires the repeated solution of eigenvalue problems with dimensions 
of the order of 10^°. At present this can only be achieved by employing elaborate numerical techniques, e.g. the kernel 
polynomial expansion method, on leading edge supercomputers.^ 

Let us first consider the adiabatic {a = 0.1) weak coupling regime Sp/t <C 1. The data presented for the single- 
particle spectral function in Fig. El give clear evidence that the system behaves like a metal. The Fermi energy 
{Kp = 7r/2) is located in the center of an only weakly renormalized band (the bandwidth is about 4t, which is the 
value for noninteracting system). The band dispersion can be derived by tracing the lowest (uppermost) peak in each 
K sector of the IPE (PE) spectra. Most notably these peaks have a spectral weight close to one. As an effect of the 
electron-phonon coupling phonon satellites separated by the bare phonon frequency wp/i occur in the vicinity of the 
tight-binding band electron levels. These excitations, however, have extremely small (electronic) spectral weight. 

If we increase the electron-phonon interaction a gap feature starts to develop in the single-particle spectra at the 
Fermi momentum. This is the situation shown in Fig. El Obviously a critical interaction strength is necessary to trigger 
the Peierls transition at finite phonon frequencies. This has to be contrasted to the result obtained for the adiabatic 
limit ujQ — 0, where a Peierls instability occurs at any finite electron-phonon coupling. FigureOalso demonstrates the 
mixed (electron-phonon) nature of the excitations. Now the spectral weight is almost equally distributed among the 
different peaks in each K sector. 

Finally we examine the behavior of the single-particle spectra in the strong electron-phonon coupling regime. As 
can be seen from Fig.^ a wide band gap emerges, indicating massive charge excitations accompanied by multi-phonon 
absorption and emission processes. The Peierls distortion of the lattice is directly connected to a charge-density-wave 
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formation (we found almost localized electrons at every second site). As a result a symmetry-broken ground state 
may occur in the infinite system, reflecting true long-range order. Now the Fermi energy is located at the center of 
the band gap and consequently the system shows insulating behavior. 

To supplement the discussion we display in Fig. [S] the contribution of states with m phonons contained in the 
ground state. From top to the bottom the coupling strength increases. At weak electron-phonon coupling Sp/t = 0.1 
the zero-phonon state is clearly the dominant one. States with more than two phonons are negligible. For the Peierls 
phase {sp/t = 1.6) the phonon distribution of the ground state follows a Poisson-like distribution, where the ratio 
Ep/'^o = 16 gives a good estimate of the mean phonon number. 

B. PRM results 

In order to integrate the renormalization equations derived in Sec. II we consider a lattice of = 1000 sites in 
one dimension. The width of the energy shell AA was taken to be somewhat smaller than the typical smallest energy 
spacing of the eigenstates of 7io,A- Fig. shows the electronic spectral functions A'^[ijj) and A~^{ijj) [calculated from 
Eqs. H35|) and (|35|l ]. and the renormalized phonon dispersion Cuq from the PRM approach for a g-M&hxe of g/t = 0.10 
[ep/t = 0.1) which is below the critical electron-phonon coupling gc- The wave number k is fixed to kp ~ 7r/2. Note 
that for g values smaller than g^ the system is in a Luttinger liquid (metallic) phase. As is seen from Fig. the 
one-particle spectrum consists of a main pole and of satellites due to the coupling of the electrons to phonons. They 
follow from the first and to the second and third terms on the r.h.s. of Eqs. (|36|l and 13 7|) and will be denoted in 
the following as coherent and incoherent excitations. Note that due to the used approximations to derive the PRM 
equations it is beyond the scope of the present approach to discuss possibly differences between a Luttinger liquid 
and a Fermi liquid behavior in the one-particle spectral functions. 

As was mentioned before we restrict ourselves to k = kp- Therefore the coherent excitation is located at the Fermi 
level [dotted line in Fig.|Bl(a)]. It is the dominant excitation of the spectrum since its pole strength is close to its 
maximum value of 1. Due to (|38|l the intensity sum over all coherent and incoherent excitations is equal to 1. The 
system is in a metallic state. Note the small finite onset of the incoherent parts in >l^(w) and A'^{uj) close to the 
Fermi level which is due to the finite phonon frequency uJq « o^o (for small 17- values). In Fig. El (b) the renormalized 
phonon frequency ujq is shown as a function of the wave number. Due to the coupling between the phononic and the 
electronic degrees of freedom, ujq has gained some dispersion. For the chosen g value this dispersion is relatively small 
except for wave numbers q close to the Brillouin-zone center and close to the Brillouin-zone boundary (fco = ±7r). This 
feature can be understood from the phonon renormalization equation H24|) . There, the second order renormalization 
contribution corresponds to a frequency dependent energy shift which is similar to that known from the phonon 
self-energy contribution of second order. The dispersion results from particle-hole excitations c^Ck+q around the 
Fermi level with a wave number transfer q which is either small (q « 0) or about q k k^ — ±tt (since at half-filling 
£kF ~ ^kp+ko)- Note that Fig.Elshows an overall agreement with the ED results of Fig.|21where the same parameter 
have been used. 

In Fig. [7| a) and b) the same quantities arc shown for g/t — 0.30 {£p/t = 0.9). In this case the coupling strength 
is somewhat below the critical value gc- This can be seen from the pole strength of the coherent pole which has 
considerably decreased from almost 1 to a value already very small compared to 1. Since the vanishing of has 
to be interpreted as a signature of the metal-insulator transition the system is still in the metallic state. In Fig. [71b) 
the renormalized phonon frequency cUq is shown. For large wave numbers approaching the Brillouin-zone boundary a 
strong softening is observed. Only for wavevectors extremely close to the boundary a sudden stiffening is found again. 
This feature is reminiscent of the behavior at ko discussed before. 

This stiffening is absent in Fig. |Hlb) where the PRM result for the renormalized phonon frequency is shown for 
g/t = 0.34 {sp/t — 1.156) which is somewhat larger than the critical g value (within the PRM approach). In this case 
the phonon excitation becomes negative very close to the Brillouin-zone boundary /cq. The negative phonon frequency 
indicates the occurrence of an instability at the critical coupling strength gc which is associated with the transition of 
the system to a Peierls state. Note that for g > gc the present PRM treatment breaks down since a possible shift of 
the ionic equilibrium positions was not taken into account. A discussion of the phonon dynamics in the Peierls state 
will be subject to further research. However, in contrast to the discussion of the phonon modes the discussion of the 
electronic excitations within the PRM method is not restricted to g values g < gc but can be extended to some small 
range above gc (see below). In Fig.jSl a) the spectral functions ^^(1^) and A^{uj) from the PRM treatment are shown 
for the same g value g/t = 0.34. The spectrum has changed from a strongly peaked distribution for smaller 3 to a 
more Poisson-like shape. The coherent pole has completely vanished af,^ — and a gap has opened in the spectrum 
at the Fermi level. The system has undergone a phase transition to a nonmetallic state. To verify that the gap does 
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not depend on the lattice size we have apphed the PRM to a rather large system with 1000 sites but also to smaller 
systems with 800, 500, 200, and 100 sites. In all cases the size of the gap was practically independent of the system 
size so that the gap can be considered as an intrinsic property of the Holstein model. As in any renormalization group 
procedure the results could depend on the actual cluster size due to additional degeneracies imposed by symmetry, 
especially for small systems. Such a dependence was not found for the systems considered here. 

To find the critical coupling strength at which the phase transition occurs we have plotted the coherent pole strength 
at wave vector k = kp for different values of g (Fig. |5Jl. Note the strong decrease of a^^ with increasing g. A 
closer inspection of the data shows that at 5 = ~ 0.31t the pole strength becomes zero which marks the transition 
of the system to the nonmetallic state. This value is somewhat larger than the critical value found from the exact 
diagonalization [gc/t « 0.24 (Sp/t ^ 0.6), compare Fig. Note, however, that in the PRM approach a rather 
large system with 1000 lattice sites was used, whereas in the ED the system had to be restricted to 8 sites. For 
the smaller systems with 800, 500, 300, 200, and 100 sites the critical coupling was slightly smaller and a pproxi mately 
g^/t « 0.30. A comparison with the critical value gc/t ~ 0.28, obtained from DMRG calculations of Refs. Tl2ll3l shows 
that the critical values from PRM might be somewhat too large. The difference is probably due to the factorization 
approximation in the PRM approach and to the simplified ansatz ()15|l for the operator form of the generator Xx.ax- 

An alternative criterion to determine the critical coupling strength is to take that value at which the renormalized 
phonon frequency Uq at the Brillouin-zone boundary q = ±ko vanishes. In this way one finds for the system with 1000 
sites gc/t « 0.30, i.e., a result which is somewhat smaller than that found from the vanishing of aj.^. Note however 
that the latter gc value seems to be the more reliable one. This can be understood from the comparison of the 
renormalization equations (|32|l for ak,\ and H24|l for ujq^x. As can be seen, for ak,\ an approximate exponential decay 
with A is found whereas a logarithmic decay follows for Wq^x- Thus, since the dominant renormalization contributions 
always occur at small A, the numerical errors are much smaller for ak,x than for Wg,A- Therefore, within the present 
renormalization approach the critical coupling strength found from the vanishing of the coherence strength oif.^ is 
probably more reliable than that given by the vanishing of (jjq=ko ■ On the other hand, from comparison with other 
approaches one has to admit that the value of gc/t « 0.31 is possibly somewhat too large probably due to the 
factorization approximation which was used to derive the renormalization equations. 

As was mentioned above, the PRM approach is only valid in the metallic regime g < gc- However, one can assure 
oneself of the fact that it may also be applied in a small parameter regime g > gc- We again consider k = kp- Due 
to the 0-functions Qk,q in all renormalization equations a renormalization approximately occurs when the energy 
difference \oJq.x + £k,\ — £k+q,\\ lies within a small energy shell between A and A — AA. Moreover, the most dominant 
renormalization processes take place for small cutoff A. Therefore, ekp+q,x ~ ^kF,\ has to be fulfilled (where a small 
phonon frequency was assumed). It follows that q « iko, where kg is the zone-boundary wave vector. According to 
(|32|l the coherent pole strength ^ will mostly be renormalized by wave numbers q « kg. In contrast /3kp,ko,x and 
lkF,ko,\ are both renormalized only once when A is small. According to akp.x becomes zero for small A in the 
insulating regime g > gc, so that the renormalization of Pkp,ko,\ and jkF,ko,>^ i^ negligible due to (|33|) and (|34|l . Note 
that the incoherent excitations in l|35|l and (|36|) close to g « fco contain negative phonon energies cuq. However, they 
do not influence the electronic spectral functions. 

V. SUMMARY 

In this paper the electron-phonon coupling in the one-dimensional Holstein model at half-filling has been investigated 
using both the projector-based renormalization method (PRM) and an refined exact diagonalization technique in 
combination with the kernel polynomial method (ED). In this system a metal-insulator transition occurs, accompanied 
by the formation of Peierls distorted state. This transition has been analyzed in terms of the (inverse) photoemission 
spectral functions, where the phonon dynamics is fully taken into account. 

Different from the ED the present PRM treatment is restricted to values g smaller or equal to the critical electron- 
phonon coupling gc- The discussion of the coupling regime g larger than gc is more complicated and has to be 
postponed to a future investigation. However, it turns out that the results for the electronic spectral functions are 
not restricted to g < gc but are also valid in a small regime regime g above gc- Therefore, the opening of an electronic 
gap can be observed. The reason is that electronic properties are rather insensitive against a small number of instable 
phonon modes close to the Brillouin-zone boundary. 

Although both the PRM and the ED technique are valid for any dimension, we have restricted ourselves to a 
one-dimensional lattice. The (inverse) photoemission spectral functions from both approaches clearly indicate a 
metal-insulator transition and the opening of a gap if the electron-phonon coupling becomes larger than a finite 
critical value. The single-particle spectra excitations are accompanied by multi-phonon absorption and emission 
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processes in the Peierls phase. The PRM provides a renormaUzed phonon dispersion which shows a softening at the 
BriUouin-zone boundary. This effect even occurs for coupling strength g smaher than the critical value so that the 
phonon softening can be understood as a precursor effect of the gap formation. 
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APPENDIX A: TRANSFORMATION OF THE OPERATORS 



In this appendix we evaluate the transformation from A to (A — AA). As an example let us consider the operator 

OO 

e^-.-4cke-^^ - = e^^ ^ (^Ck) = J] ^X';,',^;, (cj^Ck) (Al) 

n=0 

Here, a new super-operator Xa,aa was introduced which is defined by the commutator of the generator X\j\\ with 
operators A on which Xa,aa is applied, Xa.aa-4 = [ATa^aAi-^]- Thus we have to evaluate the commutators on the 
r.h.s. of HA1|I 

[^A,AA, 4ck] = Yl {^k-q:q:A0k-q,q(A, AA) (fo^jcj^.qCk + 6qcj,Ck-q) 

^ q 

-[k^(k + q)]} (A2) 

By applying X^^aa twice on c|j.Ck products composed of four fermionic and of two fermionic and two bosonic operators 
occur. In order to keep only operators which are also present in 7i(A-AA) a- factorization is performed. One obtains 

^[^A,AA,[^A,AA,4ck]] = (A3) 

= ^ II {^k-q,q,A0k-q,q(A, AA)«;, + 71^^ J 4-qCk-q 

q 

- -Bk-q,q,AQk-q,q(A, AA)(1 - Tlfc-q.A + "q^) ^Ck 

+ 5k-q,q,A0k-q,q(A,AA)K_q ;,-n^J6^5q + [k ^ k + q] } 

where the expectation values n^^ A' "q a have been defined in (|19() . The third order term of (|A1|) again gives interaction- 
type contributions. Thus, all operator terms appearing on the r.h.s. of IjAlfl are traced back to those which also appear 
in Tix- This property enables us to evaluate all higher order commutators with Xa,aa and thus the transformation 
(|A1I) . the result is given in ((TH|l . 
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FIG. 1: Low-cncrgy excitations of the single-particle spectral function A^{lli) and the corresponding integrated spectral density 
[S'lf{u)) = A'^{(jj')du)'] aX K = ±7r/2 for £p/t = 0.6 and uJo/t = 0.1. Dashed (solid) lines show the spectrum including all 
phonon modes (without the (5 = phonon mode). 



FIG. 2: (Color online) Wave-number-resolved spectral densities for photoemission [A^iyj); dashed (red) lines] and inverse 
photoemission [A^(a;); solid (black) lines] in the metallic state (sp/t <C 1). The corresponding integrated densities S^{oj) are 
given by bold lines. All ED data were obtained for an eight-site system with periodic boundary conditions. Note that the sum 
rules are fulfilled. 



FIG. 3: (Color online) PE [dashed (red) lines] and IPE [solid (black) lines] spectra near the metal to Peierls insulator transition 
point. 
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FIG. 4: (Color online) PE and IPE spectra in the Peierls phase. 



FIG. 5: (Color online) Phonon distribution in the ground state of the spinless Holstein model for three characteristic coupling 
strengths. 



FIG. 6: Results from the PRM approach for g/t — 0.10 (cp/t = 0.1): (a) Electronic spectral functions A^(a;)(black) and 
A^(lj)(t:cA) where the wave number k is fixed to kp = tt/2. The coherent excitation peak (dotted line) is at a; = 0. (b) 
Renormalized phonon dispersion Coq. The original phonon frequency is uo/t = Q.l. 
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FIG. 7: Same quantities as in Fig. |S| from the PRM approach for a g-value of g/t — 0.30 {sp/t = 0.9) somewhat below the 
critical value gc/t ~ 0.31 {sp/t — 0.96). 



FIG. 8: Same quantities as in Fig. Elfrom the PRM approach for g/t = 0.34 (sp/t — 1.16) which is somewhat larger than 
g^/t ~ 0.31 {sp/t = 0.96). Note that the coherent pole has vanished and a charge gap has opened (Fig. 8 a)) 



FIG. 9: Coherent pole strength q^^ from the PRM approach as function of the electron-phonon coupling g for a system with 
1000 lattice sites 
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